### Sample size calculation
n.for.2p <- function (p1, p2, alpha=0.05, power=.8, ratio=1) {
  if (any(p1 <1) & any(p2 <1)) {
    r1    <- ratio +1
    pbar  <- (p1+ratio*p2)/r1
    sqrt1 <- sqrt(r1 * pbar * (1-pbar))
    sqrt2 <- sqrt(ratio * (p1 * (1-p1))+ p2*(1-p2)  )
    n0    <- (((qnorm(1-alpha/2)*sqrt1) - (qnorm(1-power)*sqrt2))^2)/
      (ratio*((p2-p1)^2))
    n1    <- (n0/4)* (1+sqrt(1+2*r1/(n0*ratio*abs(p1-p2))))^2
    n1    <- trunc(n1) +1
    n2    <- trunc(ratio * n1)
  }else{
    stop("Both p1 and p2 must be less than 1")
  }
  if(length(alpha) > 1 ) {alpha1 <- alpha }else {alpha1 <- NULL}
  if(length(power) > 1 ) {power1 <- power }else {power1 <- NULL}
  if(length(ratio) > 1 ) {ratio1 <- ratio }else {ratio1 <- NULL}
  table1 <- cbind(p1, p2, n1, n2, alpha1, power1, ratio1)
  colnames(table1)[colnames(table1)=="alpha1"] <- "alpha"
  colnames(table1)[colnames(table1)=="power1"] <- "power"
  colnames(table1)[colnames(table1)=="ratio1"] <- "n2/n1"
  returns <- list(p1=p1, p2=p2, n1=n1, n2=n2, alpha=alpha, power=power, ratio=ratio, 
                  table = as.data.frame(table1))
  class(returns) <- c("n.for.2p", "list")
  returns
} 

### print n.for.2p
print.n.for.2p <- function(x, ...){
  if(nrow(x$table) < 6){
    cat("\n")
    cat("Estimation of sample size for testing Ho: p1==p2", "\n")
    cat("Assumptions:", "\n", "\n")
    cat("     alpha =", x$alpha, "\n")
    cat("     power =", x$power, "\n")
    cat("        p1 =", x$p1, "\n")
    cat("        p2 =", x$p2, "\n")
    cat("     n2/n1 =", x$ratio, "\n", "\n")
    cat("Estimated required sample size:", "\n", "\n")
    cat("        n1 =",x$n1,"\n")
    cat("        n2 =",x$n2,"\n")
    cat("   n1 + n2 =",x$n1+x$n2,"\n","\n")
  }else{
    cat("Assumptions:", "\n", "\n")
    if(length(x$alpha)==1) cat("     alpha =", x$alpha, "\n")
    if(length(x$power)==1) cat("     power =", x$power, "\n")
    if(length(x$ratio)==1) cat("     n2/n1 =", x$ratio, "\n")
    cat("\n")
    print(x$table)
  }}
